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^ i ABSTRACT 

! We investigate the atmospheres, oceans, and crusts of neutron stars accreting at rates 

sufficiently high (typically in excess of the local Eddington limit) to stabilize the burning 
of accreted hydrogen and helium. For hydrogen-rich accretion at global rates in excess 
CSJ ■ of 1O~ 8 M yr _1 (typical of a few neutron stars), we discuss the thermal state of the deep 

ocean and crust and their coupling to the neutron star core, which is heated by conduc- 
tion (from the crust) and cooled by neutrino emission. We estimate the Ohmic diffusion 
time in the hot, deep crust and find that it is noticeably shortened (to less than 10 8 yr) 
from the values characteristic of the colder crusts in slowly accreting neutron stars. As 
\ suggested by Konar & Bhattacharya, at high accretion rates the flow timescale competes 

with the Ohmic diffusion time in determining the evolution of the crust magnetic field. 
At a global accretion rate of M ~ 3 x lO~~ 9 M & yr , the Ohmic diffusion time across 
a scale height equals the flow time over a large range of densities in the outer crust. 
In the inner crust (below neutron drip), the diffusion time is always longer than the 
flow time, for sub-Eddington accretion rates. We speculate on the implications of these 
calculations for magnetic field evolution in the bright accreting X-ray sources. 
>• . We also explore the consequences of rapid compression at local accretion rates ex- 

ceeding ten times the Eddington rate. This rapid accretion heats the atmosphere/ocean 
to temperatures of order 10 9 K at relatively low densities; for stars accreting pure he- 
lium, this causes unstable ignition of the ashes (mostly carbon) resulting from stable 
helium burning. This unstable burning can recur on timescales as short as hours to 
days, and might be the cause of some flares on helium accreting pulsars, in particular 
4U 1626-67. Such rapid local accretion rates are common on accreting X-ray pulsars, 
where the magnetic field focuses the accretion flow onto a small fraction of the stellar 
area. We estimate how large such a confined "mountain" could be, and show that the 
currents needed to confine the mountain are large enough to modify, by order unity, the 
magnetic field strength at the polar cap. If the mountain's structure varies in time, the 
changing surface field could cause temporal changes in the pulse profiles and cyclotron 
line energies of accreting X-ray pulsars. 
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Introduction 



The thermal and compositional structure of the 
outer layers of an accreting neutron star is of inter- 
est to studies of accretion-induced magnetic field de- 
cay (e.g., Urpin fc Gcppcrt 1995| ), thermonuclear pro- 
cesses, and the thermal and compositional structure of 
the inner crust and core. In this paper, we consider 
objects accreting at rates high enough for stable burn- 
ing of hydrogen and helium (typically requiring super- 
Eddington rates; for a recent review see Bildsten 1998) 
and examine the thermal and compositional structure 
of the underlying ocean and crust. By "ocean and 
crust" we mean the region far below where the accreted 
matter is decelerated from its free-fall velocity and di- 
rectly underneath the hydrogen/helium burning layer 
(Fig. [j] provides a sectional cut through this part of the 
neutron star). We are therefore considering fluid layers 
that are gradually being compressed by the continuous 
addition of fluid at the top of the atmosphere. Previ- 
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ous works (Ayasli & Joss 1982; Fujimoto et al. 1984 



Miralda-Escude et al. 1990) have examined this prob- 
lem for stars accreting at globally sub-Eddington rates 
in a spherically symmetrical fashion. In this paper, we 
do not impose such a restriction and therefore param- 
eterize accretion by the local rate per unit area, rh. In 
general, neutron stars accrete at super-Eddington rates 
by means of an optically thick column; on accretion- 
powered X-ray pulsars, the local accretion rate in this 
column can easily exceed the Eddington rate by orders 
of magnitude. 

Our motivation is twofold. For hydrogen-rich accre- 
tion at global rates of order 1O -8 M0 yr _1 (the Edding- 
ton limit), we find the minimum temperature in the 
deep crust and the resulting maximum possible Ohmic 
diffusion time there. This is applicable to the brightest 
accreting objects in our galaxy (e.g. Sco X-l and the 
GX objects, which are weakly magnetic, B <C 10 12 G) 
and the Magellanic Clouds. For accreting pulsars with 
high magnetic fields (;> 10 12 G), we mostly consider the 
conditions underneath the accreting polar cap, where 
magnetic confinement of the accreted matter ensures a 
locally high accretion rate that rapidly compresses the 
ocean and heats it to temperatures of order 10 9 K at low 
densities. If this accreted matter contains an apprecia- 
ble amount of carbon (or oxygen), a violent thermonu- 
clear instability occurs and produces a flare from the 
polar cap. At super-Eddington local accretion rates, 
the recurrence time for this instability is short enough 
(hours to days) to allow observations of the flares. 
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Fig. 1. — Schematic of the outer layers of an accreting neutron 
star. Infalling matter (straight arrows) is decelerated from a ve- 
locity v ~ 0.5c at a shock and radiates a flux F GMfn/ R. Be- 
neath the shock, the accreted matter forms an atmosphere/ocean 
in hydrostatic balance where the downward flow velocity is v = 
m/p (§ pp. We indicate both the density and the physical depth 
(roughly the pressure scale height) of the underlying layers on 
the left-hand side of the diagram. The flux (wavy arrows) in this 
region is set by the steady-state hydrogen and helium burning 
{hereEnuc = X[(Q H /A) + (Q 3a /12)\ + Y (Q Za /12), see eq. |J}. 
Below the H/He burmng region, the matter forms an electron- 
degenerate ocean (§ U), where energy transport is mediated by 
electron conduction. The ocean is nearly isothermal, aiidthe flux 
is set by further compression of the matter (see eq. V.2V). At a 
density of roughly 10 9 gem -3 , the ions crystallize into a lattice. 
The thermal state of the crust and core is described in § p.3| 

We start in § || by deriving the Eulerian fluid 
equations for a plane-parallel geometry. Presuming 
steady-state hydrogen/helium burning (valid for super- 
Eddington accretion), we briefly investigate, in § ||, the 
thermal structure of the star's upper atmosphere. Be- 
cause the case of pure helium accretion is simpler, and 
because it is relevant to a number of binaries, we first 
study the steady-state burning of helium and show that 
the ashes of this burning are mostly carbon. When the 
neutron star is accreting a hydrogen-rich mixture at 
approximately the Eddington rate, the stable hydro- 
gen burning is never fast enough to consume all of the 
hydrogen before the helium ignites. Rapid proton cap- 
tures onto the rapidly proliferating seed nuclei then pro- 
duce a mix of heavy (trans-iron) elements. Although 
the final mix is quite complicated (and has yet to be 
calculated), the thermal profile is readily found. 

We next move down to the degenerate electron 
ocean, where the density is greater than 10 6 gcm~ 3 . 
The thermal flux flowing through the ocean is mostly 
generated by settling and compression within it. We 
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show in § [| how m determines the thermal profile of this 
deep ocean. Non-equilibrium nuclear reactions (mostly 
neutron emissions and pycnonuclear reactions) are an 
important heat source in the deep crust ( jHaensel & 
Zd unik 1990| ). Using a crude formulation of the crust 
heat sources, we show that most of this energy is con- 
ducted into the neutron star core, which finds a steady- 
state core temperature hot enough to balance this heat- 
ing with neutrino losses. 

In § [5], we use the thermal profiles in the crust to 
estimate the Ohmic diffusion time there. The effect of 
accretion-induced heating on the Ohmic decay of crust 



magnetic fields has previously been considered (Ro- 
mahi 199d; iGeppert & Urpin 19941; pethick & Sahrling 



19J5| lUrpin fc Geppcrt 1995| ; |Urpin fe Geppert 1996 
Konar fc Bhattacharya 1997| ) for stars accreting at sub- 
Eddington rates. In the deep crust, the local accre- 
tion rate is roughly independent of latitude; there are, 
however, a few systems that accrete globally at or near 
the Eddington rate. As a first cut, we compute the 
Ohmic diffusion time across a scale height throughout 
the outer crust (p <; 10 13 gcm -3 ) and find that it is 
greater than the time for matter to flow across a scale 
height for near-Eddington accretion rates. This implies 
that crust magnetic fields are advected to great depths 
before they are able to decay. 

Although variations with latitude in the local accre- 
tion rate may be negligible deep within the crust, they 
cannot be ignored in the upper ocean and atmosphere. 
In the case of accretion-powered pulsars, for which the 
infalling accretion stream is funneled onto a small area 
about the polar cap, the local accretion rate can be 
very super-Eddington even if the global accretion rate 
is sub-Eddington. The lack of type I X-ray bursts im- 
plies a locally super-Eddington accretion rate at the 
cap in accreting pulsars (Joss & Li 1980; Ayasli & Joss 
1982; Bildsten & Brown 1997). In § g, we show that the 
matter is tied to the field lines to very large depths and 
that spreading away from the polar cap does not occur 
until the accreted matter reaches pressures greatly ex- 
ceeding B 2 /8n (a detailed discussion is in Appendix |§]). 
The resulting high local accretion rate causes unstable 
carbon ignition at relatively low densities (§ [?]). The low 
ignition column depth and rapid accretion rate combine 
to make the resulting large "flares" from the polar caps 
recur on an hourly to daily timescale. This timescale 
differs strongly with the calculation of Woosley & Taam 
(1976), who found a roughly ten yr recurrence time for 
near-Eddington accretion rates. We discuss the observ- 
able signatures of such burning, and briefly compare 



our theoretical models to some accreting X-ray pulsars, 
including LMC X-4 and 4U 1626-67. We conclude in 
§ H with a summary of our main results. 

2. The Equations in Planar Geometry 

Because the accretion flow onto the neutron star is in 
general asymmetrical, we define a local accretion rate, 
771, as the mass accreted per unit time per unit area. 
The physics of the reactions, compression, and settling 
all depend on this quantity. The local Eddington ac- 
cretion rate is 
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Me777 M C 
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7.5 x 10 4 gcm~ 2 s~Ve 



R 



10 km 



(1) 

where ot is the Thomson scattering cross-section, c is 
the speed of light, R is the radius of the star, m u is the 
atomic mass unit, and the number density of electrons 
is pj (fi e m u ). Future references to vriEdd are for p e = 1. 
Because the ram pressure in the internal sub-sonic flow 
is negligible, the momentum equation is just that of 
hydrostatic balance, 



dP 

dz 



= ~P9, 



(2) 



allowing us to parameterize the spatial coordinate by 
column depth, y(z) = pdz', so that P(z) = gy(z). 
The heat flux flowing through the deep atmosphere (far 
beneath the photosphere, see Fig. [l]) is much less than 
the accretion flux (as GMm u /R), so radiation pressure 
is unimportant. Since accretion does not appreciably 
change the radial coordinate of isobars (see Appendix 
[X|), y is an Eulerian coordinate. We take the surface 
gravity g to be constant, since the pressure scale height, 
H = P/ pg, is much less than R everywhere in the atmo- 
sphere. Throughout this paper, we assume a neutron 
star mass and radius of 1.4M Q and 10 km, respectively; 
we also use Euclidean spacetime, with g = GM/R 2 . 
We ignore components of heat transport and pressure 
gradients parallel to the surface. Even in the case of 
accretion onto a small fraction of the surface, the atmo- 
spheric scale height is much smaller than lengthscales 
along the surface. As a result, the area of the top and 
bottom "faces" of the ocean under the accretion col- 
umn is much greater than the area along the "edges" . 
For a polar cap of radius 10 5 cm (a tenth of the stellar 
radius), the ratio of edge area to face area is ~ 10~ 2 . 

The continuous accretion of fresh material modifies 
the equations of continuity and entropy. The continuity 
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equation for a species j (with number density rij) is 



<9n, 



V- (rijv) 



„(i) 
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+ T 
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prod ' 



(3) 



Of 

where rj^ st {t^ &) is the destruction (production) rate 
of species j due to thermonuclear reactions. We pre- 
sume that all ions co-move (a good approximation at 
high accretion rates where the downward compression 
is faster than any drift time for a heavy particle). Defin- 
ing a mass fraction 
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where Aj is the baryon number of species j, and, recog- 
nizing that the accreted matter flows through constant 
pressure surfaces at a velocity v = raj p, we expand 
the continuity equation (eq. S) in (y,t) coordinates to 
obtain 



dJCj 
dt 



A jm U ( (j) 



P (J) ) 

prod / 
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In this derivation we neglect the contribution of binding 
energy to p. The entropy equation is 



ds 1 
T— = — VF 

dt p 



(0) 



where F is the flux and e is the sum of all sources and 
sinks of entropy. The flux obeys Fick's law, 



F = -KVT, 



(7) 



with K denoting the thermal conductivity. Heat 
transport is set by both radiative (primarily Thom- 
son scattering and free-free absorption) and conduc- 
tive (electron-electron and electron-ion scattering) pro- 
cesses. The total thermal conductivity is then 

K = 4aC T 3 + K (conduction), (8) 

3p/^rad 

where K ra d is the radiative opacity. Using the equations 
for hydrostatic balance (eq. J2j) and flow velocity v = 
mj p, and the thermodynamic relation 
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(9) 



we write the entropy as a function s = s(T(y, t), P(y)) 
and transform equations (|^) and (0) into the more use- 
ful forms 



dF 
dy 
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dT 
~dt 



dT 
dy 
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dT 
dy 
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V ad (10a) 
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Here cp is the specific heat at constant pressure, and 
V a d = (dl nT/d lnP) s is the adiabat. Equations (||), 



(10a), and (10b) with appropriate boundary conditions 
describe the thermal evolution of the outer layers of the 
neutron star. 



3. Hydrogen and Helium Burning in the Up- 
per Atmosphere 

We now derive the thermal and compositional state 
of the upper atmosphere for these characteristically 
high accretion rates. This determines the outer starting 
point for the derivation of the conditions deeper in the 
ocean and crust. The accreting hydrogen and helium 
typically fuse to heavier elements within an hour after 
arrival at the neutron star atmosphere. Because of /3- 
decay limitations, hydrogen burning via the CNO cycle 
is stable at rh > 10~ 2 rfiEdd and never fast enough to 
consume the accreting hydrogen before the helium ig- 
nites (|Lamb fc Lamb 1978|; [Taam fc Picklum 1979|; Fu- 



jimoto, Hanawa, & Miyaji 1981). For m > 10 riiEdd, 
helium burning occurs in a hydrogen-rich environment, 
which enhances the nuclear reaction chains and en- 
ergy release (Taam 1982, 198 5). All analyses (|Joss fc 
Li 1980|; |Ayasli fc Joss 1982| [Fushiki fc Lamb 1987b|; 



Taam, Woosley, & Lamb 1996) have shown that the he- 
lium burning in this mixed environment becomes ther- 
mally stable when the local accretion rate is, by coin- 
cidence, super-Eddington (see Bildsten 1998 for a re- 
view). It is in this stable regime that we carry out our 
calculations. 

Finding the compositional makeup of the ashes from 
steady-state hydrogen-rich accretion is complicated by 



the rp-process (jWallace & Woosley 1981 


; Champagne| 


k Wiescher 1992); |Van Wormer et al. 1994; Schatz et 


al. 1998 


), a sequence of rapid proton captures onto 



nuclei that can produce proton-rich elements beyond 
iron. The energy released is approximately known, 
however, so that we can roughly determine the upper 
atmosphere's thermal profile for steady-state burning. 
Before discussing this problem, we first consider the 
simpler case of pure helium accretion. This is relevant 
to a number of ultra-compact (P or b < 50 min) bina- 
ries (4U 1626-67, 4U 1820-30, and 4U 1916-05, Nel- 
son, Rappaport, fc Joss 1986; potentially X1850-087, 
Homer et al. 1996) and maybe the wide binary Cygnus 
X-3, if it in fact contains a neutron star. 
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3.1. Steady-State Pure Helium Burning 

There are a few neutron stars thought to be ac- 
creting nearly pure helium. For these stars, the first 
nuclear reaction is 3a — > 12 C. We use the energy 
generation rate from Fushiki & Lamb (1987a) where 
Qza = 7.274 MeV is the energy released per formed 
carbon. The radiative opacity K ra d is the sum of elec- 
tron scattering (with high temperature and degener- 
acy corrections from Paczyhski 1983) and free-free ab- 
sorption. We use the conductive opacity, K CO nd = 
3pA"(conductive)/4acT 3 , as given by Fushiki (1986). 
The degenerate electron pressure is calculated using an 
interpolation formula from Paczyhski (1983). 

The density at the location of helium burning de- 
creases as the accretion rate increases, which eventu- 
ally allows the 12 C (a, 7) 16 reaction to compete with 
the 3a reaction in consuming helium. We use Buch- 
mann's (1996) reaction rate with the fitting formula 
that yields S = 146 keV barns at 300 keV center-of- 
mass energy^. This rate is roughly 1.7 times the Caugh- 
lan & Fowler (1988) 12 C (a,j) 16 rate, and so is more 

co nsistent with nucleosynthetic constraints ( Weaver 

& Woosley 1993 ). Recent measurements ( Azuma ct 
al. 1994 ) also support this S'-factor value. In addi- 
tion to 12 C(a,7) 16 0, we also include the subsequent 
16 O(a,7) 20 Ne, 20 Ne (a, 7) 24 Mg, and 24 Mg (a, 7) 28 Si 
reactions. The reaction 28 Si(a,7) 32 S is unimportant 
for m < 40 mEdd, although for higher rates the re- 
action network would have to include it as well. For 
these reactions, we use rates from Caughlan & Fowler 
(1988) with the following adjustments: (1) following 
Schmalbrock et al. (1983), we set the uncertain factor 
between and 1 in the 20 Ne (a, 7) 24 Mg rate to 0.1; 
(2) we also set the uncertain factor between and 1 in 
the 24 Mg (a, 7) 28 Si rates to 0.1. Screening is treated 
according to Salpeter & Van Horn (1969). 

We simultaneously integrate the steady-state ver- 
sions of the continuity, flux, and entropy equations (eqs. 
||, [10b I, [ 10a |) as a function of column depth. Our 
integration starts far beneath the photosphere, where 
the outward heat flux is set by the energy released 
from nuclear reactions and the settling of the material. 
Though the flux coming through the envelope is much 
less than the accretion flux (« GMm/R), the tempera- 
ture just beneath the photosphere is still near the value 
set by the accretion flux. The flux at the top of the 



atmosphere would simply be F — {Q^ a /\2m u )rn = 
0.61 MeV(m/m u ) if all the helium burned to carbon 
and there was no energy released from gravitational 
settling. However, we are interested in the makeup of 
the products of helium burning; moreover, the heat re- 
leased from gravitational settling is non-negligible at 
high accretion rates. A hotter atmosphere has a greater 
scale height, and so at high accretion rates more energy 
is released by lowering a parcel of fluid a distance of or- 
der a scale height. To integrate the fluid equations, we 
start with a tria l va lue of the flux and integrate equa- 
tions ( 10a ) and ( lOt ) to the depth where all the helium 
is burned. At this point, we compare the flux to that 
coming from further gravitational settling of the ashes 
(see § ||) in deeper regions of the star, and iteratively 
adjust the flux. 
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1 Ther e is a typographical error in the fitting formula (Buchmann 
L99el, eq. [2] ). The terms (1 + p9 IT? /3 ) and (l+p B /T^) should 



be squared (Buchmann 1997) 



Fig. 2. — Atmospheric temperature and composition for steady- 
state accretion and burning of pure helium. The neutron star is 
of mass M = 1.4Mq and radius R = 10 6 cm, and the accretion 
rate is m = 10.0rh,Edd = 7.5 X 10 6 gcm _2 s _1 . The temperature 
is shown by the solid line in the top panel, whereas the dashed 
line displays the outward directed flux per accreted baryon. The 
bottom panel displays the mass fractions of the major elements 
present ( 12 C, ie O, 20 Ne, 24 Mg, and 28 Si) for pure helium burn- 
ing. 



When burning in steady-state, the helium is con- 
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sumed at the column depth y where its lifetime to the 
3a reaction is the time to reach that depth, t& rs y/rn. 
For the values of m considered here, burning occurs at 
a column depth y m 10 8 gcm~ 2 , which is reached in a 
time tfl < 10 3 seconds. Bildsten (1998) provides a sim- 
ple discussion of the physics of the thermal stability of 
helium burning. Linear perturbation calculations show 
that the helium burning is thermally unstable when 
T < 5 x 10 8 K in the burning region (Bildsten 1995, 
1998). From our steady-state calculations and a few 
time-dependent calculations, as in Bildsten (1995), we 
find that pure helium burning is thermally unstable un- 
til m £ 10m E dd- 

Figures ^, |3|, and |^ display the steady-state thermal 
and compositional structure of the outer layers of the 
neutron star accreting helium at rates of m/m-Edd = 10, 
20, and 40, respectively. Each figure shows the results 
for a different accretion rate. The top panel displays 
the temperature (solid line) and the flux per baryonic 
accretion rate, F/(rh/m u ) (dashed line). The flux flow- 
ing through the atmosphere increases with rh; this leads 
to increasingly hotter temperatures and lower densities 
at the location of the burning and subsequently more 
contributions from the a captures on carbon. For each 
of Figures ^, ||, and |[ the composition as a function 
of column depth is shown in the bottom panel. As we 
show in Figure |, at m = 40 m,Edd the final mixture is 
| carbon and ^ silicon. At lower accretion rates, the 
dominant nucleus made from helium burning is carbon. 
The outward-directed flux at the bottom of the atmo- 
sphere comes from the gravitational settling of matter 
beneath the helium burning region. This flux is deter- 
mined by continuing the integration of equations fllOtf ) 
and fll0b| ) to very deep parts of the neutron star, as 
described in section EJ 

3.2. Steady-State Burning of Accreting Hy- 
drogen and Helium 

At these high accretion rates, the hydrogen cannot 
burn rapidly enough to be depleted prior to helium ig- 
nition. As a result, hydrogen burning proceeds via the 
rp-process, which is accelerated by the proliferation of 
"seed" nuclei from helium burning. As a first guess, we 
presume that the hydrogen burns completely where the 
helium ignites and neglect calculating the composition 
of the ashes. These atmospheres are much hotter than 
the pure helium burning case (for the same rh) due to 
the larger flux flowing through the envelope; this leads 
to thermally s table helium burning once rh p> mEdd 
flBildstcn 1998|). 
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Fig. 3. — Same as Fig. 
20.0 rh Edd = 1.5 X 10 6 gcm" 



but for an accretion rate rh 
,-l 



We presume that the flux flowing through the atmo- 
sphere is that released from burning all of the hydrogen 
and helium to carbon, 
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(11) 



where Qh = 28MeV is the energy release per formed 
helium. This approximation yields lower limits to the 
true temperature at the location of helium burning in 
the presence of the more complex rp-process for two 
reasons. First, we presume that the reactions proceed 
only to carbon, although the residual nuclei are likely 
to be much heavier. In addition, Taam et al. (1996) 
have found (for high accretion rate type I bursts) that 
not all of the hydrogen is consumed at the location 
where the helium is depleted, in which case there is 
further energy release from electron capture on protons 
at larger depths in the atmosphere (for a recent discus- 
sion, see Bildsten fc Cumming 1998} ). Therefore, our 
thermal profiles for depths greater than where helium 
burns should be considered as lower limits. 

For accretion rates m/mEdd = 1-0, 5.0, 10.0, and 
20.0 and a helium fraction Y = 0.3 (Fig. |), we m- 
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Fig. 4. — Same as Fig. LA but for an accretion rate rh 
40.0 m Edd = 3.0 X 10 6 gcm- y s- 1 . 



tegrate equation (10b), with the flux given by equa- 
tion (|Tl]), downwards into the atmosphere until one- 
half of the incident helium has burned. At this point, 
we presume that all the hydrogen has burned as well, 
in which case the flux originating from deeper points is 
entirely caused by gravitational settling. For definite- 
ness, we took the composition to be either carbon (solid 
curve) or iron (dotted curve). The helium ignition 
depth (clearly denoted by the discontinuity in slope) is 
nearly independent of the accretion rate, yielding times 
to helium ignition of roughly 0.4 (mEdd/^) hr. At the 
highest accretion rate shown, the infalling helium ig- 
nites only 80 s after arrival onto the neutron star. The 
helium burning is stable once m J> riiEdd; the cool- 
ing rate at the temperature (5 x 10 s K) in the helium 
burning region at that accretion rate is fast enough to 
stabilize the burning. 

4. Gravitational Compression of the Ashes in 
the Deep Ocean 

We now discuss the thermal state of the degenerate 
ocean underneath the hydrogen/helium burning region. 
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Fig. 5. — Atmospheric temperature profiles for hydrogen-rich 
accretion (X = 0.7, Y = 0.3) at the super-Eddington accretion 
rates (from bottom to top) m/rraEdd = 1, 5 a 10, and 20. In the 
outer atmosphere, the flux is set by hydrogen/helium burning 
(dot-dashed curves). Where the helium is half-dcplctcd, wc join 
these curves onto the settling solutions (no burning) for an ocean 
of cither carbon (solid curves) or iron (dotted curves). The slopes 
are discontinuous at this point because the flux released from set- 
tling and compression alone is much less than the burning flux. If 
any carbon is present, it will ignite where de nuc /dT = de coo \/dT 
(see § 7.1; heavy dashed line). The iron settling solutions end at 
a lower column depth than the carbon settling solutions because 
of crystallization. At column depths greatat-than ~ 10 9 gem -2 , 
the density is approximately given by eq. jl7|). 

Upon completion of stable hydrogen/helium burning, 
the gas is hotter than 5 x 10 s K and denser than (2-5) x 
10 5 gcm~ 3 . The temperature and nuclear evolution of 
this matter under further compression depends on both 
the local accretion rate and the composition. 

First consider steady-state settling solutions in the 
absence of further nuclear burning. In this case the 
continuity equation (eq. f|) is identically satisfied, so 
our equations are 

dT F , s 

(12a) 



and 



dF 
dy 



dy 



dT 



c P - 



pK 

cpTm 



V ad 



(12b) 



dt y 

where e„ is the neutrino cooling rates from Schindcr 
et al. (1987). From these equations we define the char- 
acteristic timescale for neutrino cooling, t v = cpT/e v . 
As long as T < 10 9 K, this time is typically longer 
than either the time for heat to diffuse across a scale 
height, t t h = cpy 2 / Kp, or the time for the temperature 
to increase under adiabatic gravitational compression, 
tgh — y/riiVnd- As a result, the thermal profile of the 
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ocean is determined by the competition between i t h and 
i g h; comparing these timescales defines a local critical 
accretion rate as a function of depth, 



P K X7 - 
Wcrit = V a d 

cpy 



(13) 



The quantity m/m cr it characterizes the thermal pro- 
file of the ocean. For rh 3> m cr it, thermal diffusion is 
negligible in the time required for compression to oc- 
eo that the flow approaches adiabatic compression: 



V +- Vad- In tho opposite limit of low accretion ratoe 



(rh 77i cr it) there is plenty of time for the heat to dif- 
fuse and the flow is far from adiabatic. We can then 



rewrite equation (12b) as 



dF_ 

dy 



cpTra^ 

y 



'ad- 



(14) 



In this case V is self-consistently small. We now show 
that, in the deep ocean, m cr j t exceeds the Eddington 
limit by nearly two orders of magnitude. 

4.1. A Simple Example of the Critical Accre- 
tion Rate in the Ionic Ocean 

We first duplicate the original case discussed by Bild- 
sten & Cutler (1995). Consider a deep ocean where 
the degenerate electrons are extremely relativistic (p ;> 
10 7 g cm -3 ) and the ions (presumed for simplicity to be 
of a single species of charge Ze and mass Am u ) are not 
crystallized. The equation of state is then roughly 



P « -n e m e c 2 x + nik^T, 



(15) 



where n e (nj) is the number density of electrons 
(ions), m e is the electron mass, and x = pp/m e c = 
1.008(Zp/Al0 6 g cm -3 ) 1 / 3 is the relativity parameter 
(pf is the electron Fermi momentum). The pressure 
scale height is then 



' ' AAJ \GMm u /R 
and the density as a function of column depth is 



p w 7.6 x 10° [ 4 



y \ -3 



10 10 gcm~ 2 



cm" 3 . (17) 



The contribution from electron scattering to the con- 
ductivity is approximated by the Wiedemann-Franz 



law, 



if (electron scattering) = 



k B 2 T 



,n e v-\ (18) 



where v — v cc + v c \ is the sum of electron-electron and 
electron-ion colli sion frequenci es. We use a convenient 
fitting formulae ( Timmes 1992| ) for v cc and the formulae 
of Yakovlev & Urpin (1980) for u ei . 



The entropy per unit mass is (Landau & Lifshitz 



1980) 



^ion T" S e lec. 



A'B 



Ai 



In 



T 3/5 



m e c 2 x 



(19) 



so that V a d ~ h i and the specific heat at constant pres- 



sure (see, e.g., Cox & Giuli 1968, § 9.13) is 



cp 



' dTj „ + pT \dlnT 



P fd\uP\ 2 /<9hiP x _i 



d\np 



Am u 



I \ rn, e c z x 



(20) 



The critical accretion rate for this ocean is then 



™crit 



10 8 gcm~ 2 s _1 



5 x 10 s K 



.4 
2Z 



> ™Edd- 



(21) 



For rh <C rn cr ;t and a nearly isothermal ocean, we find 
the flux in the deep ocean by integrating equation (|l4j) , 
which yields ([Bildsten fc Cutler 1995]) 



F - 3 ™ k * T i n ( VP 



4 Am,, 



» 



Fn 



(22) 



where Fo is the flux at a fiducial column depth yo > V- 
This solution serves as a useful guide to understanding 
our more accurate solutions, as we will consider accre- 
tion rates in excess of mEdd, but well below m cr it. 

4.2. Accurate Solutions for the Compressed 
Ashes 

For the case rn-Edd < »rt <C rn cl i t , so that V <C V a d, 



we now numerically solve equations (12a) and (12b) 
for the settling of matter in the deep ocean. We in- 
clude neutrino cooling and include nonideal (ion-ion 
electronic) contributions to the equation of state. For 
the electrons' contribution to the specific heat, we use 
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Fig. 6. — Temperature (top panel) and flux (bottom panel) as 
a function of column depth for a carbon ocean formed by the 
steady state burning of accreted pure helium. The density as a 
function of column depth is approximately given by eq. (hit) . The 
accretion rates are JTi/rhsjdd = 5.0 (solid line), 10.0 (dotted, line), 
20.0 (dot-dashed line), and AO.Q-i^hort dashed line). For compar- 
ison, the flux predicted by eq. ( J22[ ) is depicted (light solid line) in 
the bottom panel for a temperature of 10 9 K. We also show the 
approximate integrated neutrino emissivity, e v y(m u /rn) (light 
dashed line), in units of MeV per nucleon for the fa = 40.0 mg^d 
solution. As is evident from this solution, when neutrino cooling 
is important, the thermal profile is sensitive to the location of 
the bottom boundary. This fact motivates our study of how the 
ocean and crust are thermally coupled to the core. 

the interpolation formula of Paczyhski (1983). The ions 
are in a liquid state (1 < V 170), where 

_ ZV/4^/\ 1/3 



1.1 



10 s K 
T 



10 8 gcm" 3 y A 1 / 3 



1/3 



(23) 



is the ion parameter for a single-component fluid. 
We use the parameterization of Farouki & Hamaguchi 
(1993) to account for the electrostatic contribution to 
the ionic energy. This parameterization, although de- 
fined only for T > 1, is, for 0.3 < T < 1.0, within 
ten percent of the considerably more cumbersome pa- 
rameterization of Hansen (1973), which is valid for 
< r <J 160. For the conditions of interest, V is always 



greater than 0.4; we therefore use the parameterization 
of Farouki & Hamaguchi (1993) for all T so that the 
specific heats are continuous. 



The numerical integration of equations (12a) and 
(12b) requires two boundary conditions. The temper- 



ature at the top of the ocean is determined by the 
solution in the hydrogen/helium burning region (§ ^|), 
whereas the flux leaving the top of the ocean (and en- 
tering the atmosphere) depends on the gravitational 
settling in deeper regions. We thus must specify ei- 
ther the temperature or flux at the bottom of the deep 
ocean. For the initial calculation, we presume that the 
flux at the bottom of the ocean is zero (but see § 4.3) 



Figure |^ depicts the resulting temperature (top panel) 
and flux (bottom panel) profiles for the ashes of pure 
helium burning under accretion rates m/rriEdd = 5.0 
(solid line), 10.0 (dotted line), 20.0 (dot-dashed line), 
and 40.0 (short dashed line). As was shown in § ^, for 
these accretion rates, the ashes will be very nearly pure 
carbon. For accretion rates below 40.0 mEdd, the flux 
roughly obeys F cx lny (eq. p^| ); this analytic solution 
is plotted (light solid line) for T = 10 9 K in the bottom 
panel of Figure |(| The slopes are steeper than that pre- 
dicted by equation (^2| ) because of ion-ion interactions, 
which increase both the specific heat and the value of 
V a d- Cooling by neutrino emission (light dashed line) 
clearly affects the highest accretion rate (40.0 rh-Edd) so- 
lution for y > 10 12 gcixr 2 and T > 10 9 K. Note that 
if the neutrino cooling term is significant in equation 
(12b), then the thermal profile is sensitive to the loca- 
tion of, and the flux at, the bottom boundary. Properly 
connecting the thermal equations of the ocean to those 
of the crust and core is the topic of the next subsection. 
We discuss, in section^], the stability of carbon burning 
in oceans this hot. 

For completeness, we also show in Figure [?] the flux 
per accreted baryon if the star is accreting a hydro- 
gen and helium mix. Because the products of H/He 
burning are uncertain in this case (see § |3|), we include 
a carbon ocean (solid lines) and an iron ocean (dot- 
ted lines), both accreting at rates rh/rfiEdd — 1.0, 5.0, 
10.0, and 20.0. For comparison, we have integrated 
both cases only to a column depth of 10 12 gcm -2 , cor- 
responding roughly to the depth where iron crystallizes. 
As a consequence, neutrino cooling does not affect the 
thermal profile, so that the flux is always positive for 
y > 10 12 gcm~ 2 . The difference in the flux for nearly 
identical thermal profiles (see Fig. ||) is due to the dif- 
ferences in specific heat between an iron and a carbon 
ocean. The ideal gas specific heat is proportional to 
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formula (Sliibazaki & Lamb 1989) 
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Fig. 7. — Flux per accreted baryon as a function of column 
depth, for accretion of hydrogen-rich matter. We show both a 
carbon ocean (solid lines) and an iron ocean (dotted lines). In 
order of increasing flux at the top (y = 10 9 gem -2 ), the curves 
are for accretion rates m/mEdd = 1.0, 5.0, 10.0, and 20.0. We 
use the boundary condition F(y = 10 12 gcm -2 ) = in both 



cases. At column depths greater_fhan 
is approximately given by eq. (|17[), 



10 gem , the density 



1/A (eq. [p0[ ). However, the stronger ion-ion interac- 
tions of the iron ocean partially offset this reduction 
in cp and also increase the value of V a d. The net re- 
sult is that for a carbon ocean, the slope of F/(m/m u ) 
is about 1.5 times greater than the slope for an iron 
ocean. 



4.3. Heating and Cooling of the Crust and 
Core 

The steady-state thermal profiles found in the previ- 
ous section assumed a zero-flux bottom boundary con- 
dition, or equivalently, an asymptotically isothermal 
temperature profile. If we were to assume a core tem- 
perature equal to this temperature, then a substantial 
neutrino flux would be generated. This alerts us to the 
facts that neutrino cooling in the crust is important 
and that neutrino emission from the core will have an 
effect on the crust and ocean even for m <; mEdd- The 



boi ndary conditions at the bottom of the ocean must 
account tor tins. 

For a core composed of superfluid neutrons and pro- 
tons, the neutrino luminosity is dominated by crust 
bremsstrahlung. For illustrative purposes, we use the 



L«4x 10 31 ergs 



Vl0 8 K/ ' 



(24) 



which demands a comparable heating of the core by 
conduction. There are two heat sources in the deep 
crust. First, the original crust of an accreting neutron 
star will be replaced in less than 10 7 yr when M > 
10~ 8 Af©yr~ 1 . As accretion pushes the crust deeper, 
the continuous compression forces the crust matter to 
endure a series of electron captures, neutron emissions, 



and pycnonuclear reactions (Haensel & Zdunik 1990). 
The heat deposited into the neutron star from these 
reactions is ps 1 MeV/nucleon and occurs at column 
depths ;> 10 14 gcm~ 2 . There is also heat generated 
directly from compression, which can at best generate 
about 0.5 MeV/nucleon. If the total electron capture 
flux, integrated over 4tt steradians, is inward-directed 
so that it diffusively heats the core, then in steady state 
this flux is equal to L u , which implies an equilibrium 
core temperature 



T r fa 5 x 10 8 K 



M 



1/6 



10~ 8 M yr" 



(25) 



However, an inverted temperature gradient is needed 
to carry the flux, so that a core in steady-state im- 
plies a temperature maximum somewhere in the ocean 
or crust. Throughout this subsection, we shall assume 
spherically symmetrical accretion. 

In order to estimate the thermal profile in the crust, 
we integrate equation (12a) for constant flux. In the 
Wiedemann- Franz approximation (eq. flcfl ), the ther- 
mal conductivity K is determined by the collision fre- 
quency of electrons with other electrons, with ions, 
and — where the ions are crystallized — with phonons 
and lattice impurities. Where the ions are crystallized 
and the temperature is greater than the Debye temper- 
ature, 



e = 1.76xl0 8 K(^f) (" / 

\ A J Vl0 10 gcm- 



1/2 



(26) 



we use the electron-phonon collision frequency obtained 



with a relaxation approximation (Yakovlev & Urpin 
1980| ), i/ph « 13e 2 k B T/(h 2 c). Although this formula 



is valid only when T J> 0, the steady-state tempera- 
ture (eq. ^5|) is so high that neglecting the T < Q 
correction does not contribute a significant error (less 
than a factor of 3) to K for p <; 10 12 gem -3 . We then 
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analytically integrate equation (12a) inward in column 
depth from yo to y to find the flux, 



F(y) 



7.0 , ]o ; ' [ ^ 

Zj 



!J0 



1/4 



(27) 

T(y)-T(y ) ergs - lcm - 2; 



i - (yo/y) 1 / 4 



where Tq and i*o are the temperature and flux, respec- 
tively, at y = y . 

A detailed solution of the thermal and compositional 
evolution of the inner crust and core of an accreting 
neutron star is beyond the scope of this paper. How- 
ever, we can use our results for the deep ocean to con- 
struct models, which will serve as a check on our choice 
of bottom boundary conditions in the integration of 
the fluid equations. We show the oceanic thermal pro- 
files from these models, for accretion rates 1.0 TOEdd and 
5.0mEdd (hydrogen-rich accretion), in Figure pL along 
with the iron settling solutions from Figure |5| (solid 
lines). 

In the first model (dot-dashed lines), we assume that 
the energy sources in the crust are concentrated at a 
single point y w . A fraction / of the generated flux is 
directed outward, while a fraction (1 — /) is directed 
inward. Because the nuclear energy released (about 
1 MeV/nucleon) is much greater than the energy re- 
leased from compression, we assume that the crust flux 
is constant away from y w and integrate equation ( |l2a| ) 
to obtain the temperature in the crust as a function of 
column depth. We used the hydrogen and helium burn- 
ing calculations described in § |3| with the assumption 
that the accreted material is processed to iron where the 
hydrogen burns in steady state. Using this assumption 
to fix the temperature at the top of the iron ocean, 
we then numerically integrate the fluid equations to 
the depth where iron crystallizes (from eqs. [^3j and 
Jl5|l this is yr ~ 5 x 10 13 gcm~ 2 ). This solution im- 
plies a flux F out flowing outwards from the crust. Using 
Font = frn(l MeV/m M ) to fix /, we then solve equation 
(27) to find the temperature at y w . From there, we 
again solve equation (^) with an inward-directed flux 
F[ n = — (1 — /)m(l MeV/m u ) to find the core tempera- 
ture T c . We then adjust / until the radiative luminosity 
entering the core balances the neutrino luminosity. Be- 
cause the thermal conductivity increases with depth, 
the flattest temperature gradient is realized by placing 
the heat source as deep as possible. We set this point 
at a column depth y w = 10 17 gcm~ 2 . As an example, 
for an Eddington accretion rate, the maximum temper- 
ature in the crust is T(y w ) = 7.0 x 10 8 K, and the core 
temperature is T c — 5.3 x 10 s K with a corresponding 



neutrino luminosity L„ = 8.8 x 10 35 ergs -1 . About 2% 
of the heat generated in the crust diffuses upward into 
the ocean; although not inconsequential, the change in 
temperature does not significantly change the results 
of our integration of the settling solutions. Including 
the contributions of Umklapp and impurity scattering 
to the conductivity only steepens the temperature gra- 
dient between y w and the core, so that T(y w ) would be 
even hotter. 

In the second and third models (which are some- 
what unrealistic), we suppose the crust to have no heat 
sources other than compression. Because of the high 
conductivity, V is much less than V a d throughout the 
crust (see § 4.1); nevertheless, the core neutrino lu- 
minosity demands an inward-directed flux. To esti- 
mate the effect of compressional heating in the crust on 
the equilibrium core temperature and the oceanic ther- 
mal profiles, we performed two int egrat ions, one with 



a constant flux (integration of eq. 12a |; dotted lines) 



and one with a constant temperature (integration of 
eq. [l4| ; dashed lines) crust. The constant-flux crust 
is completely devoid of any heat sources, while in the 
isothermal crust model, heat is generated via compres- 
sion. For both cases, we numerically integrate over the 
ocean, and then integrate analytically either equation 
(12a) or equation ( [Til ) to estimate the heat flux into, 
and the temperature of, the core. We again balance the 
flux directed into the core with the neutrino luminos- 
ity, and adjust the flux at the top of the ocean until 
the two match. For a constant-flux integration at m = 
"T-Edd, the core temperature and neutrino luminosity 
areT c = 3.0 x 10 s K and L v = 2.3 x 10 34 erg s"\ respec- 
tively. The isothermal crust model produces a slightly 
higher temperature and luminosity: T c = 3.3 x 10 8 K 
and L y = 4.8 x 10 34 ergs -1 . 

As a check of our model assumptions, we compute 
V /V a d at the top of the crust (where this ratio will be 
greatest). Because of the strong sensitivity of L v to 
temperature, the core temperature is roughly the same 
in all cases. We also compute the ratio of the settling 
flux (from eq. |l4}|) generated in the crust to that en- 
tering the core (in the model where the flux is assumed 
constant, we calculate the settling flux that would have 
been generated). For sub-Eddington accretion rates, 
either approximation will serve — the gradient is very 
sub-adiabatic and the settling flux in the crust is small 
compared to that required by the neutrino luminosity. 
For higher accretion rates, V/V a d is no longer negligi- 
ble, so that the constant-flux approximation is better, 
as most of the settling flux is generated in the upper 
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of heat into the core (as well as different core mod- 
els), but again, their typical crust temperatures were 
less than 10 8 K. Like Ayasli & Joss (1982), their mod- 
els had heat sources only in the outer envelope. Non- 
equilibrium reactions in the crust were first incorpo- 
rated into accreting neutron star models by Miralda- 
Escude et al. (1990). They did not consider local ac- 
cretion rates greater than roughly 0.01 times the Ed- 
dington rate (using their radii), so that the crust tem- 
peratures were typically never greater than 10 s K. 

In conclusion, the physical processes deep in the 
crust do not appreciably affect the thermal profiles in 
the upper ocean, which is dominated by compressional 
heating and nuclear processing of hydrogen and helium. 
It is important therefore to treat the upper atmosphere 
correctly, as the thermal profiles of the ocean depend 
completely on the physics in that region. 



Fig. 8. — Thermal profil es in the ocean for various models 
of the crust/core. The star is accreting hydrogen-rich (X = 
0.7, Y = 0.3) matter at rates m/»riEdd = 1-0 an d 5.0. For 
each accretion rate, we compute three different models: (1) a 
constant-flux crust with heat generation from electron captures 
(dot-dashed lines); (2) an isothermal crust with no heat sources 
other than compression (dotted lines); and (3) a constant-flux 
crust with no heat sources (dashed lines). For comparison, we 
also display the settling solutions from Figure B (solid lines), for 
which F(y = 10 12 gcm — 2 ) = 0. At column depths greater than 
~ 10 9 gcm~ 2 , the density is approximately given by eq. ([tv])- 



ocean (see eq. Q). This regime is unlikely to be im- 
portant, as the accreted matter will have spread around 
the stellar surface at this depth, so that the local ac- 
cretion rate is not more than a few times Eddington. 



Previous authors (Ayasli & Joss 1982; Fujimoto et 



al. 1984; Miralda-Escude et al. 1990) also considered 



the total energy balance of the neutron star core and 
its influence on the outer crust, albeit for much lower 
accretion rates. Because of the rapid compression of 
the star's outer layers and the reactions in the deep 
crust, the outer layers of the neutron star do not ap- 
proach an isothermal temperature profile. In contrast, 
Ayasli & Joss (1982) considered densities greater than 
10 s g cm~ 3 to be part of an isothermal core. Their cal- 
culation of neutrino luminosities differs from ours in 
that they include general relativistic corrections to the 
core temperature (their standard model was for a star 
of mass 1.41 M© and radius 6.57 km, so that the rela- 
tivistic corrections are quite nonnegligible). The cor- 
rections imply a higher neutrino luminosity, as mea- 
sured by an observer at infinity, for a given core tem- 
perature. Fujimoto et al. (1984) considered the flow 



5. Ohmic Diffusion in the Deep Ocean and 
Crust 

Present uncertainties in the composition of matter 
after hydrogen/helium burning prohibit a calculation 
of the subsequent chemical evolution of the ocean for 
hydrogen-rich accretion. However, even though the 
composition is not well known, we can still use the 
thermal profiles as estimates of the crust temperatures. 
This is important to the evolution of the magnetic field, 
as the accretion-induced heating of the crust reduces its 
conductivity and hastens the Ohmic diffusion of crust 
magnetic fields ( Geppert fc Urpin 1994 ). This heating 
also increases the mass of the ocean. These effects have 



been considered 


(Romani 1990: 


Geppert & Urpin 1994; 


Pethick & Sahrling 1995; 


Urpin & Geppert 1995; 


Urpin 


k Geppert 1996; 


Konar & Bhattacharya 1997) for stars 



accreting at M <; lO" y M yr _1 . 

There are, however, a few neutron stars accreting 
globally at or near the Eddington rate. There are 
two X-ray pulsars (LMC X-4 and SMC X-l) and the 
six bright "Z" sources (Sco X-l, GX 5-1, GX 349+2, 
GX 17+2, GX 340+0, Cyg X-2). The accreted material 
will have spread over the surfaces of these star for col- 
umn densities > 10 14 g cm -2 (see § ||), so that a spheri- 
cally symmetrical approach is warranted for this calcu- 
lation. We thus use our solutions for the thermal profile 
of the deep crust at accretion rates rh ~ ?riEdd (see sec- 
tion 4.3) to estimate the Ohmic diffusion timescales in 



the deep crust of these neutron stars. 
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5.1. The Microphysics in the Crust 

The conductivity in the crust is set by electron- 
phonon and electron-impurity scattering. In the 
relaxation-time approximation, the conductivity is 



(Yakovlev fc Urpin 1980) 



e 2 n R 



,\J\ + x 2 v 1 



(28) 



where v is the sum of the electron-phonon (f p h) and 
electron-impurity (z/; mp ) collision frequencies, for which 
we use the expressions from Urpin & Yakovlev (1980), 



1 

^ P h 
1 



h 2 c 



13e 2 /c B T 

3ttH 3 Z 1 
8m e e 4 Q x 



9 
3T5T 



-1 1/2 



(29a) 



(29b) 



Here Q = nj~ x ^AZj — Z) 2 nj parameterizes the im- 
purities, described by charge Zj and number density 
nj. We assume that x ^> 1 and neglect anisotropies 
in the relaxation time due to the magnetic field. From 
these conductivities we then calculate the local Ohmic 
diffusion time over a scale height, 



Tdiff 



A H2 

c z 



(30) 



We are using the pressure scale height H (eq. p6[ ) 
as the characteristic lengthscale. At neutron drip, 
H/R « 0.01(2Z/A) 4 / 3 (p/10 11 gcm- 3 ) 1 / 3 , so that a 
plane-parallel approach is valid throughout the crust. 

5.2. Ohmic Diffusion Times in the Crust 

For the temperature s in the crust, we used the es- 
timated profiles from § O for the case of nonequilib- 
rium nuclear reactions occurring deep in the crust. As 
in that section, we assume that temperature is only a 
function of depth y, as at these depths the accreted 
matter will have spread around the star. We plot in 
Figure || (Q = 1.0) the local Ohmic diffusion time 
(solid lines) for accretion rates of 0.5, 1.0, and 5.0 times 
Eddington. We also show the flow time over a scale 
height, tfi = y/m (dashed lines). A few conclusions 
are immediate. First, where the ions vibrate classically 
(T p> O), the ratio Tdie/tn is nearly independent of 
depth until near neutron drip; moreover, for accretion 
rates rh/m Edd < 0.23(A/2Z) 2 (T/5 x 10 8 K), the diffu- 
sion time is always greater than the time for matter to 



flow through one scale height. Second, impurity scat- 
tering is unimportant throughout the crust for Q <J 1. 
Because we placed the heat sources at a fixed depth 
U = Uw the thermal gradient changes sign there (see 
§ ^J^). Electron captures remove pressure support and 
therefore decrease the pressure scale height, causing the 
abrupt decrease in the Ohmic timescale (solid line, Fig- 
ure ^) . Once neutron pressure dominates the equation 
of state, the scale height again increases with depth. In 
this region, the flow timescale is always longer than the 
diffusion timescale for m < ?nEdd- 
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Fig. 9. — Ohmic diffusion in the crust as a function of column 
depth for accretion rates m/rh,Edd = 0.5, 1.0, and 5.0. The con- 
ductivity includes both electron-phonon scattering and electron- 
impurity (Q = 1.0) scattering. We show the timescale for Ohmic 
diffusion (solid line) over a scale height and the timescale for 
the crust to be pushed through a scale height, y/m (dashed 
line). The two timescales become comparable, above neutron 
drip, when M fs 3x 10 — 9 Mq yr — 1 . Below neutron drip, the flow 
timescale is always longer than the diffusion timescale for sub- 
Eddington accretion rates. The density, as a function of column 
depth is approximately given by eq. (p. 7|) . 

Our calculation has several differences with previ- 
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ous calculations. First, the flow timescale is always 
less than the Ohmic diffusion timescale, so that the 
downward advection of current is nonnegligible for 
M > 3 x lO _9 M yr -1 . This agrees with the find- 
ings of Konar & Bhattacharya (1997), who found that 
the decay of the crust field decreased as the accretion 
rate increased (they considered accretion rates of less 
than 10" 9 Af Q yr" 1 ). Second, the impurity scattering 
does not make a significant contribution to the con- 
ductivity until near-nuclear densities are reached (for 
Q = 1.0). Of course, there is no physical reason why 
Q should have a value near unity; in fact, it is not 
clear that the crust can be adequately described as 
a crystalline lattice with impurities. Given the non- 
equilibrium nuclear processes in the ocean and crust, a 
better description may be that of an alloy. The chem- 
ical profile from Haensel & Zdunik (1990), which used 
M = 10~ 10 Moyr -1 , assumed a single species present 
at each column depth, starting with a uniform iron 
ocean. This is unlikely to be realized in practice, es- 
pecially if the ocean is composed of rp-process ashes. 
However, increasing the number of species present in 
the crust is unlikely to raise the conductivity, so that 
our calculation probably represents an upper bound to 
the actual diffusion time. A final difference is that our 
oceans are much deeper — the melt surface is located at 
i/~5x 10 13 gcm~ 2 (p ~ 10 10 gem -3 ). 

In computing the diffusion and flow timescale, we 
use the local value of the pressure scale height H as 
the characteristic lengthscale. In contrast, Geppert & 
Urpin (1994) used a fixed lengthscale of 832 cm (which 
is much less than the pressure scale height in the deep 
crust, eq. |lq|). Because t^s is proportional to the 
square of the lengthscale, the diffusion timescale with 
their length is much shorter than ours at depths where 
the local scale height exceeds their fiducial value. 

The way in which these calculations are applied to 
the question of global magnetic field decay depends on 
where one places the currents that provide the mag- 
netic dipole. Geppert & Urpin (1994) and others have 
placed all of the current in the crust and argued that 
the surface field will decay on the timescale given by the 
Ohmic diffusion time at the depth where the majority 
of the current is flowing. Given an initial distribution of 
currents in the crust, the evolution of the surface field 
is then determined by the competition between the de- 
cay rate at the depth of the current and the downward 
advection of the currents to higher conductivity regions 
(Konar & Bhattacharya 1996). Estimates made in this 
way are very sensitive to the initial current distribution, 



as is clear from the change in the local Ohmic diffusion 
time with depth (see Fig. |)J). If, however, the field 
penetrates the core, then the diffusion time is roughly 
R/H ~ 10-100 times longer than in the absence of field 
penetration ( Pethick fc Sahrling 1995| ). 



6. 



Magnetic Confinement of Accreted Matter 
at a Polar Cap 



The most common sites for locally super-Eddington 
accretion rates are the magnetic polar caps of accreting 
pulsars. The local rates can be quite large and most 
estimates of the polar cap area lead tom> ?7iEdd at 
the photosphere for the majority of bright accreting 



pulsars (Lamb et al. 1972; 


Arons & Lea 1976; 


Eisner 


k Lamb 1977 




Bildsten & Brown 1997 


). We now find 



the depth at which the matter begins to spread away 
from the magnetic polar cap. Regions underneath this 
spreading depth have a local accretion rate closer to the 
average, M/AirR 2 . We start by showing that the flow 
on the polar cap occurs under nearly MHD conditions 
for local rates m ;> O.lmEdd- 

6.1. The Magnetic Reynolds Number in the 
Atmosphere and Ocean 

The magnetic pressure is always much larger than 
the matter pressure until rather large depths are 
reached. We thus begin in the upper atmosphere, 
where B 2 ^> 8irP, and ask whether the flow can cross 
the field lines via Ohmic diffusion. The fluid behav- 
ior is parameterized by the magnetic Reynolds number 
lZ m = ifl/Vdiff, where Tdiff is the Ohmic diffusion time 
across a scale height (defined in eq. ) . In the absence 
of lateral flow, the flow time scale over a pressure scale 
height is ifl = y/rh. In the upper atmosphere, above 
and near the hydrogen/helium burning , we estimate 
the conductivity as (see, e.g., |Spitzer 1962j ), 



2(2fc B r) 3/2 

7r 3/2 m l/ 2 e 2yY 



(31) 



where A ss 8 is the Coulomb logarithm^ at typical at- 
mospheric temperatures. We then find that the mag- 
netic Reynolds number in the upper atmosphere is 



(32) 



/ m \ / 10 5 g cm 3 
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VmEdd/ 



T 



10 s K 



5/2 



A 



2 We follow the convention of Yakovlev & Urpin (1980) and denote 
the Coulomb logarithm with A instead of In A 
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When the hydrogen/helium burni ng is in steady-state, 
we can use the flux equation (eq. 12a] ) with Thomson 
scattering opacity to eliminate p in favor of T and to 
and obtain 



I atmosphere 

1400 



(33) 



\ m Edd/ V 



1 x 5 x 10 8 K\ 1/2 /8 



A 



T 



where we have used the fiducial temperature for the 
hydrogen/helium burning region. This clearly shows 
that, for to > 3 x 10 _2 TOEdd, Ohmic diffusion is neg- 
ligible prior to stable hydrogen/helium burning. For 
lower accretion rates, the matter will spread via Ohmic 
diffusion prior to hydrogen/helium ignition. 

Moving from the upper atmosphere to deeper re- 
gions, we now consider an ocean composed of H/Hc 
burning ashes. The conductivity here is mostly set by 
electron-ion scattering, which we treated in § |^. For 
this calculation, we replace A i by unity. We are in- 
terested in the case where the ocean is composed of 
light elements, such as carbon, which might later ignite 
unstably. We then find that the magnetic Reynolds 
number in the deep ocean is 



(34) 
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m-Edd 



and is nearly independent of the density. As above, 
this serves to define a critical accretion rate, to > 
3 x 10 _3 TOEdd (where 1Z m ~ 1), above which we can 
presume nearly ideal MHD. At rh = TOEdd, the density 
at which carbon ignites is roughly 10 8 gcm~ 3 , which 
implies that TZ m w 54. The magnetic flux is therefore 
frozen into the fluid down to depths where carbon ig- 
nition occurs for the accretion rates of interest to us. 
Hence, if the magnetic field is sufficiently strong, the 
fluid will remain on the polar cap until reaching these 
depths. We now estimate how strong the field must be. 

6.2. A Simple Magnetostatic Mountain 

Polar cap accretion causes transverse pressure gra- 
dients that are balanced by curvature of the magnetic 
field. We give an order-of-magnitude estimate of the 
bending required to support a mountain of fixed over- 
pressure (Appendix [b| contains a more detailed discus- 
sion of this problem.) In this sense, we are envision- 
ing accumulation of matter, followed by eventual recon- 
nection or Ohmic dissipation. Our crude criterion for 



where spreading occurs is where the field lines are bent 
by a large angle, as originally discussed by Hameury et 
al. (1983). 

For simplicity, consider a poloidal, azimuthally sym- 
metrical magnetic field, B = (B CT ,0, B z ), at the po- 
lar cap. The overpressure due to the accumulated ac- 
cretion flow distorts the field from B = (0,0, Bq) to 
the perturbed configuration B = (f? ro ,0, B — 8B Z ). 
Because the flow velocity for m ~ mEdd is so small 
[v « 0.1 cms _1 (10 8 g cm _2 /y)], we neglect terms con- 
taining the velocity in the time-independent momen- 
tum equation, 



- WP + pg + i(Vx B) x B = 0. 

47T 



(35) 



We discuss the solution of equation (p5|), under some 
simplifying assumptions, in Appendix B. Here, we re- 
produce a simple estimate ( Hameury et al. 1983| ) of the 
value of Bcj/Bq. The characteristic length in the in- 
direction is the scale height H and the characteristic 
length in the tn-direction is the polar cap radius i? ca p- 
The current density induced by the transverse pressure 
gradients in the mountain is 



4tt 



-Vx B 



(36) 



and is toroidal. The radial component of equation (35) 
implies that 



Bo i B 



R 



cap 



4tt V H 
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cap 



(37) 



The components SB Z and B^ are related by V- B = 0, 
so B^jH w SB z R cap and equation (p7|) becomes 



AttHP 



H 



R 



cap 



H 



2R, 



(38) 



cap 



where f3 = 8ttP/B 2 is the ratio of the gas and mag- 
netic pressures. If H m 10 3 cm and R, 



cap 



10 b cm, 

then equation ( |38p implies that the field lines become 
appreciably distorted where B^/Bq ~ 1, or 



R c 



H 



100. 



(39) 



Because the magnetic field lines are bent where j3 ~ 
100, we presume that spreading occurs via recon- 
nection or interchange instabilities, and that lateral 
spreading of the matter occurs at a pressure P J> 
(B 2 1 '871") (i? cap / H) . Using the relativistic degenerate 
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equation of state for the pressure (eq. JL5|| ) and scale 
height (eq. fig] ) in equation (|38|), we find that the field 
lines will be appreciably distorted at a column 



2/sprcad 
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i\ 4/5 



(40) 
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As we show in the next section, stars accreting pure 
helium at m ^ lOrh-Edd ignite carbon unstably before 
the material can spread. 

The induced current density needed to distort the 
field creates at the polar cap a dipole moment, 



H = i /" r'xJ(r')d 3 r'. 



(41) 



From equations (J36J) and (|39|), we see that this dipole 
moment opposes the external field and scales as 



H 



2i? c 



BnR, 



O-ttcap 

~8^ 



(42) 



Though the polar cap is radially very thin, the induced 
dipole distorts the magnetic field topology for a dis- 
tance i? cap above the surface and therefore can eas- 
ily modify the flow of matter there. A totally self- 
consistent calculation that includes the effects of this 
induced moment on the accretion flow is beyond the 
scope of this paper, but it is interesting to note that 
mountains at these large overpressures (if stable) can 
easily change the local magnetic field strength in the 
polar cap region. Whether this can cause observable 
temporal changes to the pulse profile or cyclotron line 
energies is not presently known. In Appendix ^, we 
describe our numerical solution to equation (|35) that 
confirms the estimates of the spreading criterion (eq. 
]39|]) and the induced dipole moment (eq. |42|). 

7. Unstable Thermonuclear Ignition of Carbon 

Carbon is the predominant element left from the 
steady-state (or even unstable) burning of pure helium. 
Taam & Picklum (1978) first showed that pure carbon 
oceans never burn in steady-state but rather undergo 
unstable thermonuclear flashes when pycnonuclear ig- 
nition occurs at densities w (6-9) x 10 9 gcm~ 3 . For 
local accretion rates near Eddington, this would result 
in a recurrent explosion every ten yr or so, as originally 
envisaged by Woosley & Taam (1976). We now em- 
phasize that, for high accretion rates, unstable carbon 
ignition can occur at much lower densities, leading to 



more frequent deflagrations. We first discuss the car- 
bon burning thermal instability and then apply it to a 
few accreting objects. 



7.1. Thermally Unstable Carbon Burning 



For the 12 C + 12 C reaction, we use the thermonuclear 
reaction rate, (av), from Ca ughlan fc Fowler (1988) 



with a strong screening factor (Ogata, Ichimaru, & Van 
Horn 1992 ). The depth where the energy loss from neu- 
trino cooling equals the energy generation from carbon 
burning is shown in Figure [w| (light dashed curve); to 
the right of this curve heating from carbon burning is 
faster than neutrino cooling. Carbon burning is still 
thermally stable above this curve, however, since ther- 
mal conduction is so efficient. Hence, unlike in collaps- 
ing white dwarfs and massive stellar cores, the ther- 
mal stability of carbon ignition depends on the com- 
petition between nuclear heating and thermal diffu- 
sion. We therefore define the unstable carbon igni- 
tion curve (Fig. [HJ heavy dashed line) by the equation 
(de nuc /dT) P = {de cool /dT) P ( |Fujimoto ct al. 1981 ; 
Fushiki & Lamb 1987b), where e coo i = pKT/y z is a 
local representation of the cooling rate from diffusion. 

Upon reaching the unstable ignition curve, a fluid 
element containing carbon will heat up on a nuclear 
burning timescale. A preliminary investigation of the 
thermal stability (both in the linear regime and a brief 
series of non-linear calculations under the approxima- 
tions outlined in Bildsten 1995) determines that the 
thermal instability ensues soon after carbon burning 
initiates. The cooling rate from thermal diffusion [con- 
tours of constant thermal time (Fig. |l^, light dotted 
lines) are shown for log(tth/ls) = 2.0, 3.0, and 4.0] 
is faster than the heating rate prior to ignition. Once 
ignition occurs, the entire thermal profile moves up- 
wards until the heating rate greatly exceeds the cool- 
ing rate. Timmes & Woosley (1992) showed that pure 
carbon burning fronts move at speeds greatly in excess 



of 10 km 



and have intrinsic widths much less than 



1 cm. As Bildsten (1995) noted, combustion fronts with 
widths much less than a scale height propagate very 
easily and are difficult to quench, as all the nuclear en- 
ergy released goes into lateral (not vertical) heating. 
Thus, a local instability should lead to ignition of all 
connected and combustible carbon layers in the star 
and a large energy release. 
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Fig. 10. — Temperature profiles (solid curves) for the case of 
pure helium accretion. We take the ashes of helium burning to 
be pure carbon. In order of increasing temperature, the solutions 

Siown areTor rh/rh^dd = 5.0, 10.0, 20.0, and 40.0 (see Figs. bl 
, kL and H). Carbon ignites where de nnc /dT = de coo \/dT {heavy 
dashed curve). The filled circles indicate where nuclear burning 
has decreased the abundance of helium to half of its initial value. 
To the right of the curve e u = e n uc (light dashed line), carbon 
burning heats the ocean faster than neutrino emission can cool 
it. We show contours of constant thermal diffusion time (light 
dotted curves) for log(t t h/ls) = 2.0, 3.0, and 4.0. At column 



depths grcatcii_Lhan 
given by eq. (|17j> . 



10 9 gcm 2 , the density is approximately 



7.2. Recurrence Times and Energetics 

As shown above, a pure carbon ocean under com- 
pression at high local accretion rates initiates ther- 
monuclear burning prior to electron capture reactions, 
which occur at p = 3.9 x 10 10 gcm~ 3 . We use the 



asymptotically isothermal solutions developed in § 4.2 
to find the ignition conditions as a function of the local 
accretion rate. This treatment is adequate, as the ther- 
mal structure of the upper ocean is mostly determined 
by the steady-state helium burning and compression, 
not by processes in the deep crust and core. In order 
of increasing temperature, the solutions (Fig. [l(], light 
solid lines) are for rh/ifiEdd — 5.0, 10.0, 20.0, and 40.0 
(see Figs. [|, ||, ||| and ^|). The filled circles indicate 
where nuclear burning has decreased the abundance of 
helium to half of its initial value. It is the intersec- 
tions of these curves with the carbon ignition curve 
(heavy dashed line) that sets the characteristic recur- 
rence time, and hence the energy released during the 
instability. 

For example, the solution for rh/rriEdd = 5.0 in- 
tersects the ignition curve at y — 1.85 x 10 11 gcm -2 
(see Fig. [l(]). The recurrence time is therefore 5 x 
10 5 s, and the burst energy (assuming that fuel ac- 
cumulated over AttR 2 of the surface burns) is 1.7 x 
10 41 (Af/10~ 8 M o yr~ 1 )erg. The thermal time, which 
sets the duration of the resulting burst in the ab- 
sence of convection, at this depth and temperature 
(7.1 x 10 8 K) is i t h = 8.4 x 10 3 s. In contrast, the solu- 



tion for fn/rh-Edd = 40.0 intersects the ignition curve at 
(y = 3.8 x 10 9 gcm- 2 , T = 1.1 x 10 9 K), so the recur- 
rence time is 1.3 x 10 3 s, implying an energy release of 
4.6 x 10 38 (Af/10~ 8 Ai© yr^ 1 ) erg; the thermal time here 
is t t h = 470 s. The energy released in the instability is 
roughly the energy released per gram of reactants times 
the amount of material accumulated between bursts. It 
is possible that not all of the accumulated matter will 
be burned at once; for example, if the accretion geom- 
etry is two disjoint polar caps, then each cap will burst 
independently. What is accurately known is the time- 
averaged ratio of the energy released in quiescence to 
the energy released during a burst (the a-parametcr). 
For carbon burning, this ratio is a ~ 105. 

7.3. Comparison to Observations 

On what neutron stars can this instability occur? In 
low magnetic field objects, the matter will have spread 
around the star prior to carbon ignition, in which case 
all of the weakly magnetic stars accreting pure helium 
(4U 1820-30, 4U 1916-05, and X1850-087) will have 
to < TO-Edd at the depth of carbon ignition. These stars 
will then undergo ignition at such high column den- 
sities (> 10 13 gcm~ 2 ) that the time to reach ignition 
will be ~ 10-100 yr, resulting in released energies of 
order 10 44 erg(M/lO~ 8 M yr-^^/lO yr). No such 
events have ever been identified. This is not surprising, 
given the long recurrence time and short duration of 
the event. 

Another possible binary with helium-rich accretion 
at high accretion rates is Cygnus X-3. van Kerkwijk 
et al. (1996) have recently tightened the case that this 
X-ray binary consi sts of a WN Wolf-Rayet star orbiting 
a compact object ( |van Kerkwijk et al. 1992 ; van Kerk- 
wijk 1993 ). Neither pulsations nor type I X-ray bursts 
have ever been seen from this object. Detailed X-ray 



timing and spectr a l studies with EXOSAT (Berger fc 



yan dcr Klis 1995 ; Smit fc van der Klis 1996 ) did not 



allow for a simple classification of this object in terms 
of the standard accreting neutron star phenomenology 
(i.e., Z and atoll sources). We would expect the helium 
burning to be unstable if the neutron star is accreting 
nearly pure helium at a rate less than ten times the 
Eddington limit. It might be possible that the type 
I bursts are obscured by the well-documented corona 
about this object or that the compact object is a black 
hole. The carbon burning would be unstable on this 
object and potentially have shorter recurrence times of 
order years. 

The most intriguing scenario for our work is an ac- 
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crcting X-ray pulsar, where the high magnetic field con- 
fines the accreted matter onto a small fraction of the 
stellar surface. The enhancement of the local accre- 
tion rate then decreases the recurrence time, possibly 
to values as short as a day, even at globally low accre- 
tion rates. The methodology by which unstable nuclear 



builmng on an accreting X-ray pulsar is used to place 
constraints on trie surtace magnetic held and polar cap 
size is discussed for GRO J 1744-28 in Bildsten & Brown 
(1997). 

For one pulsar, 4U 1626-67, our calculations are di- 
rectly applicable. Chakrabarty et al. (1997) inferred a 
global accretion rate in excess of 2 x 10 _10 M Q yr _1 , 
and the constraints on the orbital parameters (P rbit — 
42 min) i mply that the accreti ng companion is a helium- 
rich star ( Chakrabarty 1998|) . This pulsar rotates at 



-Pspin = 7.68 s and has a corotation radius of 6.5 x 
10 s cm. If we presume that the magnetospheric radius 
is near equilibrium, so that the accretion follows the 



field lines from the cor otation radius ( |Lamb ct al. 1973 
Ghosh fc Lamb 1979] ), then the local accretion rate 



at the polar cap (each polar cap has area 10 cm ) 
is 10 6 gem -2 s _1 . Carbon ignition then occurs at a col- 
umn depth yi gn ss 1.9 x 10 10 gem -2 , which the accreted 
material reaches within 4 hr of deposition onto the sur- 
face. The ignition occurs long before the matter has 
spread around the star for B ;> 8 x 10 11 G. The ther- 
mal time at that depth is roughly 100 s. These quanti- 
ties set the recurrence and duration of unstable carbon 
burning for this particular polar cap model. Unfortu- 
nately, the recurrence time is very sensitive to the local 
accretion rate (and hence is a sensitive function of the 
accretion geometry and global rate). In this case it is 
better to invert the problem: given a recurrence time, 
we can immediately infer the local accretion rate, and 
thus the polar cap area. From Figure [LOL the recurrence 



time is roughly 6 x 10 4 (10.0mEdd/^) s. There have 
been many reports of flares on shorter timescales from 
this pulsar (recurrence times of 1000 seconds). These 
flares could result from carbon burning if the matter is 
either confined to an area smaller (roughly 2 x 10 9 cm 2 ) 
than our naive estimate, or if the global accretion rate 
is a factor of ten higher than the minimum stated ear- 
lier, as suggested by the recent distance estimate of 
Chakrabarty (1998). If the polar cap is larger than 
10 10 cm 2 , as suggested by accretion models in which the 
matter penetrates the magnetopause via a Rayleigh- 
Taylor instability (Axons & Lea 197C; Eisner & Lamh 
1970), then these flares could be caused by unstable he- 
lium burning. The recurrence time would then imply a 



local accretion rate of order 10 5 g cm -2 s -1 (low enough 
for helium burning to be unstable), or a polar cap area 
10 11 cm 2 , which is still much less than the total 



4 

surface area. 



We were initially very interested in the large daily 
flares from the X-ray pulsar LMC X-4 (Levine et 
al. 1991|; Dal Fiume et al. 1997). During these ~ 30 min 



flares, the pulsed fraction increases, the pulse becomes 
single peaked and sharply defined, and the spectrum 
softens. The ratio of the persistent fluence to the 
time-averaged burst fluence (the a parameter) is ~ 33; 
hence, the energetics of the flares are reasonably consis- 
tent with the unstable burning of accumulated matter. 
The time-averaged luminosity of rs 4 x 10 38 ergs -1 im- 
plies a super-Eddington global accretion rate of approx- 
imately 4 x lO _8 M0yr -1 , and so we investigated the 
possibility that the flares could be from carbon burn- 
ing. The difficulty is that the massive companion is 
providing hydrogen-rich material, in which case we do 
not expect large amounts of carbon in the deep ocean, 
as the hydrogen burning is known to produce much 
heavier ashes. It presently seems unlikely that carbon 
burning is the cause of these flares. 

Even though we cannot presently explain them, 
we would like to note three other accreting pulsars 
that exhibit large flares: SMC X-l, GX 304-1, and 
EXO 2030+375. EXOSAT observations of SMC X-l 
(Angelini, Stella, & White 1991) show an - 80s flare 
with a rise time ~ 1 s. This behavior is qualitatively 



duplicated in GX 304-1 QMcClintock ct al. 1977| ), which 
has flares of duration ~ 100 s and rise time ~ 10 s. The 
flares of EXO 2030+375 ( farmer et al. 1989| ) are char- 
acterized by both a longer rise time (~ 30 min) and a 
longer duration (~ 1-2 hr). Moreover, the recurrence 
time is ~ 4 hr. Because EXO 2030+375 is a transient X- 
ray pulsar, it is expected that the accretion rate varies 
greatly in time. With the launch of RXTE, it is pos- 
sible to obtain detailed pulse and burst profiles from 
these sources. For example, the sharpening of the pulse 
profile on LMC X-4 during a burst would be expected 
from a single polar-cap burst. If this is true, then it 
should be possible to determine if only one polar cap 
bursts, or if the polar caps alternate. These types of 
observations would help in the eventual unraveling of 
the cause of these flares. 



8. Summary and Conclusions 

We have investigated the thermal and compositional 
structure of a rapidly accreting neutron stars in many 
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different environments, both highly and weakly mag- 
netic. We first demonstrated, in § 3.1, that for accre- 



tion rates rh J> lO.OmEdd (f° r which pure helium burn- 
ing is stable), the ashes of pure helium burning are 
mostly carbon (Figs, g, [| and ||). In § |7.1|, we calcu- 



lated the unstable carbon ignition curve; the intersec- 
tion of our settling solutions with that curve (see Fig. 
flO| ) determine the recurrence times and energetics of 
the unstable carbon burning. For the high local accre- 
tion rates underneath a pulsar polar cap, we found that 
episodes of unstable burning can recur within hours to 
days. This is much shorter than the previously found 
recurrence times of years for lower accretion rates. In 



7.3 , we listed some promising sources for observing 
such unstable carbon burning. We are especially in- 
trigued by the 1000 s flares on 4U 1626-67, which ac- 
cretes from a helium-rich star. There have also been 
large flares seen from LMC X-4, SMC X-l, GX 304- 
1, and EXO 2030+375. Further observations of the 
burst energetics, the recurrence time, the burst dura- 
tion, and the changes in the spectrum, the pulse profile 
and phase, and the pulsed fraction of the luminosity 
during a burst should reveal the nature of these events. 

Accreting X-ray pulsars offer a promising laboratory 
for studying the effects of a locally high accretion rate. 
In the upper atmosphere and ocean, the Ohmic diffu- 
sion time is always much longer than the flow timescale 
(§ |6.l| ), so that the matter is tied to the field lines. We 
showed, in § |6.2| , that the magnetic field keeps matter 
from spreading (which occurs when the initially vertical 
magnetic field is bent so that the poloidal component 
is -B ro ~ Bo) around the star until the gas pressure is 
P ~ (Rcap/H)B 2 /8ir. Where the gas pressure is suffi- 
cient to bend the field lines, a dipole moment of order 
-Bo-Rcap 3 is induced (see Appendix |b|). This induced 
dipole is in opposition to the stellar dipole and can af- 
fect the magnetic field strength, and consequently the 
accretion stream, at the polar cap. If time dependent, 
this could lead to observable changes in the pulse pro- 
files or cyclotron line energies from these objects. 

We also disc ussed the case of hydrogen and helium 
burning (§ 3.2) for accretion rates where the burning 
is stable. We described in detail the thermal state of 
the ocean, crust, and core for neutron stars accreting 
globally at rates in excess of ~ 10 _9 A/ Q yr _1 . In § 
we used these thermal profiles to estimate the Ohmic 
diffusion time throughout the star. We found that the 
hotter crust of a rapidly accreting neutron star short- 
ens the diffusion time to less than 10 s yr throughout 
the crust for an Eddington accretion rate. Moreover, 



for accretion rates rh ~ mEdd, the flow time is less than 
the diffusion time throughout the crust. This raises in- 
triguing questions about the evolution of global mag- 
netic fields in rapidly accreting neutron stars, like the 
bright "Z" sources. 
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A. Transformation to Column Depth Coordi- 
nates 

There is a subtlety in the plane-parallel equations for 
subsonic accretion flow: namely, whether constant pres- 
sure contours (isobars) are spatially fixed or move radi- 
ally outward during accretion. (For a brief discussion, 
see § 3.2 of Bildstcn, Salpeter, & Wasserman 1992.) For 
example, accretion of an incompressible liquid would 
lead to outward motion of the isobars (relative to the 
stellar center of mass) at a rate Vi SO = rh/ p. In this ap- 
pendix we show that allowing for compression implies 
that the isobars have constant altitudes on timescales 
shorter than that needed for the accreted matter to 
reach pressures where the equation of state becomes 
nearly incompressible. For neutron stars this timescale 
is greater than 10 6 yr, i.e., the equation of state is very 
compressible until the internuclear spacing is on the 
order of a nuclear diameter. 

The transformation equations between plane- 
parallel and column depth coordinates are 
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is the upward rate of motion of an 
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isobar i>i S0 . For example, consider a polytropic gas of in- 
dex n, P cx sitting on an incompressible floor at 
z = 0, and let quantities evaluated at z = be denoted 
by a subscript 0. Hydrostatic balance, dP/dz — —pg, 
implies (we are neglecting the ram pressure) 



P = Po 



1 



(A2) 



where A = gpo/[(n + 1)Pq] is the total height of the 
polytropic atmosphere. The column depth y — y(z, t) 
is 



y = A/?o- 
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(A3) 



Now, Po/g equals the accumulated column mass of 
the atmosphere, yo- Accretion will clearly increase yo 
(dyo/dt — m). Differentiating equation ( |A3| ) with re- 
spect to t while holding y fixed yields 



Po 

so that the velocity of a fluid element is 



(A4) 



in I p 



(A5) 
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v = h v iso 

P P V Po 

The advective derivative is then simply 

— —-— ' (l-£\ — 
dt dz dt \ po J dy 

d . d . . „. 

~8t +m d-y (A6) 

far above the incompressible region where p <C po- 



B. Magnetic Confinement of the Accreted 
Matter: A More Detailed Approach 

To investigate the magnetic confinement of a mound 
of accreted matter in more detail, we considered the 
following idealized problem. In a cylindrical geometry 
(vu, </>, z), we parameterize the magnetic field by a func- 
tion ^(b, z) such that 



B = —Vijj x ba,, 
w 



(Bl) 



where is the basis vector in the azimuthal direction 
and ip(0, z) — 0. We are assuming that the magnetic 
field is poloidal, i.e., B = (B m ,0, B z ). This parame- 
terization of B ensures that V- B — is satisfied iden- 
tically and that the flux <!> crossing a disk of radius R 



in a plane of constant z is $ = 2m/)(R, Zq). As a re- 
sult, field lines are tangent to lines of constant ip. This 
technique is well-known in solar physics, where this ap- 
proach (with some additional assumptions about the 
nature of the toroidal field) leads to the Grad-Shafranov 
equation (see, e.g., Freidberg 1987), a variant of which 
we now derive. 

Hameury et al. (1983) performed a similar calcula- 
tion and concluded that no significant bending of the 
field lines occurred. This approximation is only true if 
f3 (H/R cap ) 2 <C I, where H is the height of the region. 
In their numerical calculations, Hameury et al. (1983) 
used integration heights of 30 and 70 m; the correspond- 
ing values of [(3(H / L) 2 ] z=0 were » 0.7 and » 60. Thus 
their approximations were not valid at the bottom of 
the region of integration, where the greatest deforma- 
tion of the field lines occurs. 



Using the parameterization in equation (Bl), we 
transform Vx B into 



Vx B 



d 1 dtp 
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In this problem, we are considering the gas to be 



polytropic, P = Up i - 1+1 / n ^ 
w = J dP/p = (n + 1)P 
the pressure in equation 
rameterization of B (cq. IB 



pV(w + V) 




so that the enthalpy is 
Using this to eliminate 
and substituting the pa- 
we have 
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VipA^p = 0, 



(B3) 



where V is the Newtonian potential. The component of 
equation (BS ) parallel to a field line is F(ip) = w(ip, z) + 
V(z) = constant. Inverting this equation then gives 
the density as a function of ip and z, 

1 



p(ip,z) 



(n+l)H 



TO) - v{z)] 



(B4) 



Thus, F(ip) specifies how each field line crossing the 
boundary of the cylinder is loaded with matter. The 
component of equation (B3) perpendicular to a field 
line implies that 

dF 1 

dip 



(A^) = 0. 



(B5) 



The system of equations (B4) and (B5) is then closed 
by specifying F(tp) and boundary conditions [at w = 
(0, -Reap) and z — (0, h)] on tp. We choose tp so that, 
along the boundary, the magnetic field matches a uni- 
form field. 



(B6) 
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In other words, all of the field lines enter and leave the 
computational box and are presumed to match onto a 
uniform field above and below the region of integration. 
The choice of F(ifj) is somewhat arbitrary in the ab- 
sence of detailed information about the accretion flow. 
A convenient parameterization is 



Polar Cap Radius [cm] 
5x10" 10 6 



F = Fn 



exp 



B m 



cap 



(B7) 



with this choice we can simulate either a hollow or filled 
polar cap geometry. In addition, we require that the 
solutions have a topology where all field lines exit the 
boundaries of the region. This ensures that F(ip) re- 
flects the field configuration for matter arriving from 
z — > ao. 

As a test of the crude scaling argument (eq. p8|), 
we solved equation (B5) for a region of one scale 

height with F = F [exp(-V>/£o-Rca P ) + !]• The to P 
of the region was chosen to have fixed density, p = 
2 x 10 6 gcm -3 . We varied independently the scale 
height (by changing the surface gravity), the radius of 
the polar cap, and the magnetic field strength. We 
chose values for which the distortion of the field was 
slight. The value of B^/Bq scaled linearly with h and 
L, and inversely as Bq~ 2 , in agreement with equation 



(38) 



In Figures Bll and B12 we show two numerical so- 
lutions to equation (p35|). Panel a depicts the constant 
density contours (thin lines) and the lines of constant 
flux (thick lines). In panels b and c we show the inte- 
gral of J over the vertical coordinate z and the radial 
coordinate w, respectively. The matter overpressure in- 
duces a dipole moment, opposed to the stellar field, of 
order B R cap 3 inside the region. The existence of large 
current sheets at the top and bottom of the region is 
evident in panel c. These sheets are required by the 
boundary conditions that the field be uniform outside 
the region. In Figure Bll (Fig. B12), the parameter 
/3(F/i? cap ) varies from 0.24 (0.97) at the top to 0.37 
(6.1) at the bottom. As can be seen in Figure B12 a, the 
maximum displacement of a field line from its equilib- 
rium location is of order 500 cm, so that B^/B z is of or- 
der unity. (The field lines are bent more than is visually 
apparent from panel a due to the different scales for the 
vertical and lateral displacements). The induced dipole 
(not including the current sheets at the top and bot- 
tom) is fi = -0.078B o i?ca P 3 e, (fi = -0.82B R cap 3 e z ) 
for Figure |Bll| (Fig. |B1~2| ). 

While the field lines threading the crust are anchored 
in place, the field lines through the outer atmosphere 




5x10" 10 5 

Polar Cap Radius [cm] 

Fig. Bll. — The confinement of accreted matter by a magnetic 
field. In panel a we show isochores (thin lines) for a matter 
distribution, at z — > oo, of F oc [1 + exp(— |ro/R cap | 2 )], where 
w is the coordinate along the surface. The magnetic field (heavy 
lines) is assumed to be uniform in the region outside the plot. In 
panels b and c we show the current density integrated over z and 
w, respectively. The parameter /3(H / -Reap) varies from 0.24 at 
the top to 0.37 at the bottom, so that no appreciable spreading 
occurs, even though /3 = 8itP/B 2 = 74 at the top (z = 0). The 
induced dipole moment (not including the current sheets at the 
top and bottom) is /x = — 0.078-BoRcap 3 e z . 

are free to move outward. Because of the large induced 
dipole, the field lines in the atmosphere will be dis- 
torted from the matter-free configuration for a distance 
~ -Reap 3> H above the ocean. 
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